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ERROR ANALYSIS OF A WEIGHTED LEAST SQUARES PROCESS TO 
DETERMINE THE LUNAR GRAVITATIONAL FIELD 

By Robert H. Tolson and Harold R. Compton 
Langley Research Center 

SUMMARY 

An error analysis has been made to determine how accurately the coefficients of 
the spherical harmonics of the lunar gravitational potential function can be estimated by 
using earth-based range and range-rate measiarements of a lunar satellite. Parametric 
analyses were made to study the effects of certain of the orbital elements on the standard 
deviations of the gravitational coefficients. The covariance matrices in which the stand- 
ard deviations are given were obtained by forming and inverting the normal matrices 
in a weighted least squares process. Standard deviations of some of the gravitational 
coefficients were strongly dependent on the semimajor axis, inclination, nodal position, 
and eccentricity of the lunar satellite orbit. Also some of the high correlations associ- 
ated with the estimates of the coefficients were found to be dependent either on the length 
of the tracking-data arc or on certain of the orbital elements. The two data types, range 
and range-rate, having weights of 10 meters and 0.002 meter per second, were found to 
be similar or essentially the same data types in the sense that they produce similar cor- 
relation matrices. There were no significant differences in the values of the determinants 
of these matrices, and combining the two data types did not lead to significant reductions 
in high correlations. From the analyses it appears that the best tracking schedule would 
be one in which the allowable tracking time is distributed over various times of the month. 
When tracking data from lunar satellites having different orbital inclinations are com- 
bined to estimate the gravity coefficients, the condition of the normal matrix for inversion 
is improved and high correlations between the even zonal coefficients are reduced. 

INTRODUCTION 

An important contribution to the scientific community would be an improvement in 
the present knowledge of the lunar gravitational field. This information could be used in 
the planning of future manned missions to the moon. For instance, close lunar satellites 
will be used to take pictures of the moon which, in turn, will aid in the choosii^ of landing 
sites for the manned spacecraft. An accurate knowledge of the lunar gravitational field 
is necessary for calculating a precise ephemeris of the photographic satellite which is 



needed for correct interpretation of the photographic data. Precise information about 
the lunar gravitational field would also aid those concerned with the figure of the moon 
and its internal structure. 

The gravitational field of the earth has been approximated through analysis of the 
trackir^ data from close earth satellites and it is expected that similar analyses of 
tracking data from close lunar satellites can be used to estimate the lunar gravity field. 
An approximation of this gravity field can be made by determinir^ a sufficient number of 
coefficients, Cn^m expansion of the lunar gravitational potential func- 

tion U in spherical harmonics: 


U =^| 
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where fx is the gravitational constant of the moon, r is the distance from the center 
of the moon to the satellite, (p is the selenographic latitude of the satellite, X is the 
selenographic longitude of the satellite, R is the mean radius of the moon, and Pj^ 
are the associated Legendre polynomials (n is degree and m is order). The coeffi- 
cient Cq q represents a deviation of the mass of the moon from a nominal value assumed 
in p. 

At the present time, knowledge of the lunar gravitational field is very limited. A 
brief discussion is given in reference 1 concerning present estimates of the mass of the 
moon and the two second-degree gravity coefficients, C 2 q ^2,2* Uncertainties in 

the current estimates of C 2 q and C 2 2 nre given as 1 x 10“^ or 2 X 10“^ in a nominal 
value of -2.071 x 10“^ for C 2 g ^^nd about 1 x 10“^ in a nominal value of 2.072 x 10“® 
for C 2 2* Estimates, such as those given in references 2 and 3, have also been made of 
some of the higher degree coefficients. It appears, however, that present estimates of 
the gravity coefficients are not very accurate and it is of interest to analyze a method of 
determining the coefficients to a higher degree of accuracy. 

The purpose of this paper is, therefore, to perform an error analysis in which a 
method for determining the gravity coefficients is analyzed. Two general approaches to 
the problem of determining the gravitational field were considered for use in this error 
analysis. These approaches are referred to as the short period or direct method and the 
long period and secular method. Both methods may be used in a differential correction 
process. In the direct method, the partial derivatives of the observables with respect to 
the coefficients are obtained from a direct formvilation of the observables as functions 
of the coefficients. In the long period and secular method, the observables or tracking 
data types are used to determine sets of intermediate elements. These elements are then 
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used as observables in a differential correction process where the long period and secular 
variations in the elements are analyzed to determine the gravity coefficients. See ref- 
erences 1 and 4 for further discussion of the distinction between these two methods. 

The gravitational coefficients are usually referred to as either zonal or tesseral 
coefficients. The zonal coefficients are those with m = 0 and the tesserals are those 
coefficients with m 0. A subset of the tesseral coefficients, where n ** m, are some- 
times referred to as sectorial coefficients. The zonal coefficients in the gravitational 
field of the earth are best determined by the long period and secular method, whereas a 
direct method is the only way of determining the tesseral coefficients. In the determina- 
tion of the gravitational field of the moon, however, both the tesseral and zonal coefficients 
are determinable by the long period and secular method. (See ref. 1.) Therefore, for a 
determination of the Ivmar gravitational field, a choice is available between the direct 
method and the long period and secular method. In this paper the direct method of deter- 
mining the gravity coefficients is analyzed. 

A parametric error analysis was made to investigate the correlations between the 
estimated gravity coefficients, the condition of the normal matrices for inversion, and 
the accuracies of estimating the gravity coefficients. The coefficients were never actually 
determined and in particular only the covariance matrix for the coefficients was needed. 
This covariance matrix was obtained by forming and inverting the normal matrices in a 
weighted least squares process. 


SYMBOLS 

A matrix containing partial derivatives of a given data type with respect to the 

gravitational coefficients 

a semimajor axis of lunar satellite orbit 

Cn rn gravitational coefficient (n is degree of the coefficient and m is the order) 

D distance from center of earth to center of moon 

DET determinant of correlation matrix 

e eccentricity of lunar satellite orbit 

i inclination of orbital plane of lunar satellite to earth-moon plane 
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m order of a particular gravitational coefficient 

n degree of a particular gravitational coefficient 

Pn associated Legendre polynomials (n is degree, m is order) 

q order of correlation matrix 

R mean radius of moon 

r distance from center of moon to lunar satellite 

®n m gravitational coefficient (n is the degree of the coefficient and m is the 

order.) 

t^ time of periapsis passage 

U gravitational potential function defined in equation (1) 

V true anomaly of the lunar satellite 

W weighting matrix 

X,Y,Z coordinates axes with origin at center of earth (The X-axis is positive in the 
direction from the center of the earth to the center of the moon, the Y-axis 
is positive in the direction of rotation of the moon, and the Z-axis is posi- 
tive in such a direction that it forms a right-handed axis system.) 

parameter to be estimated (a particular gravity coefficient); i = 1,2, . . . j 

Aa covariance matrix of estimated parameters 

X selenographic longitude of the lunar satellite 

Amin minimum eigenvalue of correlation matrix 

Ajj^^ maximum eigenvalue of correlation matrix 

p gravitational constant of moon 
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f2 longitude of ascending node of lunar satellite orbital plane measured in earth- 

moon plane in direction of rotation of moon from positive X-axis 

longitude of ascending node of lunar satellite orbital plane at the beginning of 
the tracking interval 

o) argument of periapsis, angle measured in lunar satellite plane from ascending 

node to periapsis 

4> selenographic latitude of the lunar satellite 

p range or distance from center of earth to position of lunar satellite 

p range rate or radial velocity of lunar satellite with respect to center of earth 

correlation coefficient; denotes correlation between aj and a 2 

a standard deviation or one-sigma uncertainty (When this symbol appears with 

a subscript, it is taken to mean the one-sigma uncertainty of the variable 
indicated by the subscript.) 

standard deviation or one -sigma uncertainty in an observation or data t5T)e 

6 = 0 ) + V 


GENERAL DISCUSSION 

The tracking geometry for a close lunar satellite is considerably different from 
that of a close earth satellite. For instance, close earth satellites pass over the tracking 
station very rapidly, and, consequently, the station can view the satellite for only a limited 
time and over only a short segment of the orbit. In general, for orbits with high inclina- 
tion, an equatorial tracking station has at most two viewing periods per day as the rotation 
of the earth carries the station through the orbital plane. Hence in order to obtain a good 
sampling of the gravity field and a suitable determination of the gravitational coefficients, 
tracking data must be obtained from a number of tracking stations located at various lati- 
tudes. For a close lunar satellite, an earth -based tracking station can view a l\mar satel- 
tile from moon rise to moon set except for the time when the satellite is occulted by the 
moon. Thus, with exception for the time the satellite is occulted by the moon, two tracking 
stations properly situated on the earth can obtain almost complete coverage of each orbital 
period of the l\mar satellite. 
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Even with such complete coverage, one might anticipate difficulty in trying to deter- 
mine coefficients with tracking data which cover a short time span because in this case 
the tracking geometry is approximated by the so-called Stationary- Moon geometry. The 
Stationary-Moon geometry is defined as the case when the line from the tracking station 
to the center of the moon remains inertially fixed during the tracking interval. It has 
been shown that with this tracking geometry the state of the lunar satellite cannot be 
completely determined. (See ref. 5.) In reality the moon is not stationary and the line 
from the tracking station to the center of the moon is not fixed but rotates about the 
Z-axis at a rate of about 0.55 degree per hour. Therefore in practice, due to this small 
angular rate, the tracking geometry approaches that of the stationary moon when short 
data arcs are considered. 

In reference 6 it was found that, after 1 orbit of tracking, the normal matrix formed 
to solve for the state of the lunar satellite could not be inverted with single-precision, 
8-decimal arithmetic; however, double -precision, 16-decimal arithmetic, proved to be 
adequate. Hence one would not be surprised to encounter similar difficulties in trying to 
determine the gravity coefficients from short tracking intervals, and it is of interest to 
investigate the length of the tracking interval required for a first determination of the 
coefficients. This will be discussed in a subsequent section of this paper. 

For a close lunar satellite, the distance from the tracking stations to the satellite 
is extremely large compared with that of a close earth satellite, 380 000 kilometers as 
-compared with 4 200 kilometers. The hour angle and declination of the satellite can be 
measiored by earth-based radar to about 0.183 degree. (See ref. 7.) These measure- 
ments lead to a very poor determination of the position of the satellite and these angular 
measurements cannot be used effectively in either the orbit-determination process or the 
process of determining the gravity coefficients. Therefore, only range and range-rate 
measurements were considered for use in the present analysis. 

ANALYSIS 

Certain simplifying assumptions were made in order to reduce the complexity of 
the equations involved. First, it was assumed that the moon revolved about the earth in a 
circular orbit. The earth-moon plane was chosen as the fundamental plane of reference 
and all the orientation angles of the lunar satellite orbital plane are given relative to this 
fundamental plane. The geometry of the problem is illustrated in figure 1. Throughout 
this analysis only one tracking station was considered; this station was taken to be at the 
origin of the geocentric coordinate system shown in figure 1. Range and range -rate 
measurements were simulated from this station with the assumptions that all measure- 
ment errors of a given data type were uncorrelated, unbiased, and of equal weight. Here 
"measurements were simulated" means that the measurements were only assumed to 
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have been made. No real tracking data were used nor were any tracking data generated 
because, in this type of analysis, neither is necessary, as will be shown subsequently. 
Standard deviations of range and range rate were assumed to be independent of time and 
to have values of 10 meters and 0.002 meter per second, respectively. These accuracies 
are estimates of those applicable to the NASA deep space network (DSN) tracking system. 
(See ref. 8.) Note that, by assuming the tracking station to be located at the geocenter, 
it is never occulted by the earth; however, occultation of the satellite by the moon was 
accounted for. In reality, a topocentric tracking station would sometimes be occulted by 
the earth and tracking would be lost as the station went over the horizon, but, in general, 
another station would come up and it would be possible to continue tracking without loss 
of coverage due to occulation of the tracking station. This is the same type of coverage 
which was assumed by having the trackiiig station at the geocenter. In general, a given 
topocentric tracking station would be able to view a lunar satellite orbit at slightly dif- 
ferent angles because of the parallax of the earth and this would normally result in a 
better estimate of the orbit and the gravity coefficients. Therefore, the standard devia- 
tions of the estimates of the gravity coefficients presented herein are considered to be 
conservative because of the assumption of a geocentric tracking station. 

The gravitational coefficients were assumed to have been estimated by using the 
simulated range and range-rate measurements in a weighted least squares process. Nor- 
mal matrices were formed and inverted in this process to obtain the covariance matrix 
associated with the estimates of the coefficients. Because the measurements were 
assumed to be uncorrelated and of equal weight, the weighting matrix was diagonal with 
equivalent terms on the diagonal, and hence the covariance matrix for one data type can 
be written as 


A„=(aTwa)-1.o2JaTa)-‘ = 


4 . 
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I^Q!1Q!2 Q?1 0^2 


Pa. at. '^a. ^a. 
1 ] 1 ] 


Pa^a^W 

o2 

“2 


0 (JO. 
Q!2 Q!j 


Pa^afa^^ 




( 2 ) 


7 



I III II III I ■■■lll■l■■ I I null I III ■ I 


where (a*^WA) is called the normal matrix, A is the matrix containing the partial deriv- 
atives of the data type with respect to the coefficients to be determined, and is the 

standard deviation of the data t 3 rpe. The standard deviation of the estimates of the coeffi- 
cients is denoted by cr subscripted with the appropriate variable whereas the correla- 
tion coefficients are denoted by p subscripted with the appropriate variables. When the 
covariance matrix for the simultaneous use of the two data t 5 ^es is desired, one simply 
adds the normal matrix for each data type and then inverts the sum. No actual measure- 
ments are ever needed in this process; only the statistics of the measurements and the 
partial derivatives in the A matrix are required. Because of the special form of the 
weighting matrix, it can be seen from equation (2) that, for one data t 3 q)e, the standard 
deviations of the estimates of the coefficients are proportional to the standard deviation 
of the data type. 

In subsequent discussions the term "correlation matrix" will be referred to fre- 
quently. This is a symmetric matrix having diagonal elements equal to 1 and the correla- 
tions between the gravity coefficients as the off diagonal elements. The elements of this 
matrix are dimensionless. 

For the direct method referred to previously, the observations are related directly 
to the gravitational coefficients given in equation (1) through the integrals of the equations 
of motion of the spacecraft. Hence, in order to form the partial derivatives of the meas- 
urements with respect to the gravitational coefficients, which are needed to form the 
normal matrix, the equations of motion of the lunar satellite must be integrated. For 
this part of the analysis, a first-order general-perturbation method similar to that dis- 
cussed in reference 9 was utilized to obtain the desired integrals of motion. The only 
difference between the equations of motion in reference 9 and those used in this analysis 
was in the choice of the element used to represent the position of the satellite in the orbi- 
tal plane. The perturbed mean anomaly is used in reference 9 whereas the mean anomaly 
at epoch is used in this analysis. The necessary partial derivatives were then obtained 
by the so-called chain rule. First, the derivatives of the observables with respect to the 
state of the satellite at each observation time are obtained from purely geometrical rela- 
tions. These derivatives are given in reference 6. Next, the derivatives of the state with 
respect to the gravitational coefficients at the observation times are obtained by direct 
differentiation of the integrals of motion. The partial derivatives of the observables with 
respect to the coefficients of interest can then be obtained by applying the chain rule to 
the appropriate derivatives. 

To make a parametric study of the effects of a given element on the accuracy of 
determining the gravitational coefficients, five elements of a nominal orbit were held con- 
stant and the sixth was varied over a given range. One exception to this procedure was 
the eccentricity variation, where, instead of using the nominal value for the semimajor 
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axis, a value of 5 000 kilometers was used. A large semimajor axis was necessary to 
insure that, over the range of eccentricities investigated, the pericentron of the orbit was 
never less than the radius of the moon. Nominal values of the set of elements used in this 
analysis are: 

a = 2 686 kilometers 
i = 15° 

Oq = 30° 
w = 0° 
e = 0.336 
tjj = 0 second 

In this study, the effects of varying w and tg on the accuracy of determining the 
gravitational coefficients were not considered. It is recognized that a variation in oj 
would cause a different portion of the lunar satellite orbit to be occulted by the moon from 
the tracking station. However, it was found from occultation studies that, for orbits with 
small to medium eccentricity, the covariance matrices did not change significantly with 
changes in a>. Also, whenever the orbital element o> is used in the partial derivatives, 
it always appears as an angle added to the true anomaly in the argument of either a sine 
or cosine function. Because the true anomaly rotates through 360° each orbit, the value 
of the sine or cosine changes through 1 period, regardless of the value of w. Hence, the 
effects of a variation in w on the accuracy of estimating the gravitational coefficients 
are negligible as long as integral orbits of tracking are considered. Likewise, the effects 
of a variation in the orbital element t^j are negligible for integral orbits of tracking 
because the parameter is used only to determine the location of the satellite at any given 
time in the orbit, and hence the measurements which are equally spaced in time over the 
orbit are independent of tg. 


RESULTS AND DISCUSSION 

All results presented herein are for the simultaneous estimation of a set of 11 grav- 
itational coefficients. This set includes the mass of the moon, all the second-degree 
coefficients, the zonal coefficients through degree five, and two fourth-degree second- 
order tesseral coefficients. Unless otherwise noted normal matrices were formed by 
assuming that 26 equally spaced observations, less the points omitted because of occulta- 
tion, were made each orbit over a Ivmar satellite orbit having the nominal orbital elements 
given in the previous section. The 26 data points were usually reduced by 15 to 20 per- 
cent because of occultation of the satellite by the moon. As mentioned previously, the 
standard deviations in the range and range-rate measurements were assumed to be 
10 meters and 0.002 meter per second, respectively. 


9 



Effects of Tracking Time 


The variation of the standard deviation of estimating the set of gravitational coeffi- 
cients with the number of orbits tracked is shown in figure 2. Here the normal matrices 
for each data type and for the simultaneous use of both data types have been summed and 
inverted after each orbit of tracking. After about 3 orbits of tracking, it appears that 
range-rate data have a slight advantage over range data as far as the accuracy of deter- 
mining the coefficients is concerned. However there are no very large differences 
between the results for either range or range rate. The simultaneous use of range and 
range rate reduced the standard deviations by a factor of approximately ^l2, as might be 
expected, since the number of data points was doubled. Of the 11 coefficients estimated, 
the two odd zonals were determined with the best accuracy because these two coefficients 
have no high correlations with the other coefficients considered. The standard deviations 
of the remaining 9 coefficients decrease very rapidly because of the reduction in the cor- 
relations between certain coefficients as the tracking interval increases. In general, 
lower correlations usually result in lower standard deviations and a better separation of 
the coefficients in the solution vector. 

The correlation matrices associated with the estimation discussed in the previous 
paragraph are shown in figure 3 after the first and tenth orbits. The high correlations, 
that is correlations greater than 0.90, have been divided into two categories. First, there 
are those correlations due to the short data arc, and, second, those correlations which 
are strongly dependent on certain of the orbital elements. In the first category, the cor- 
relations are very high, 0.99 or larger, after 1 orbit of trackii^ and decrease to reason- 
able values after a few orbits of tracking. Several examples of this type of correlation 
can be found in figure 3. For instance, in the range-rate matrix the correlation between 
C 2 1 and S 2 j is -0.9992 after only one orbit and then decreases to -0.42 after the 
tenth orbit. In the second category there are correlations which are strongly dependent 
on the nodal position and inclination and to a lesser extent on the eccentricity. These high 
correlations do not change much over a tracking time of a day or so. The nodal-dependent 
correlations are usually between pairs of tesseral coefficients whereas the inclination- 
dependent correlations are usually between the zonal coefficients. Examples of these 
correlations are pointed out in subsequent discussion. 

One of the primary concerns in any real operational estimation procedure is the 
condition of the normal matrix for inversion. The determinant of the correlation matrix, 
which must be greater than 0 and equal to or less than 1, gives some indication of the 
condition of the normal matrix. (See theorem 4.1, p. 14 of ref. 10.) It has been a common 
practice to analyze the condition of the covariance matrix on the basis of the ratio of max- 
imum to minimum eigenvalues. If this ratio is equal to 1, the matrix is diagonal and 
hence can be inverted trivially without loss of significant figures. However, the eigen- 
values of the covariance matrix are dependent on the dimensionality of the parameters 
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being estimated, and therefore a diagonal covariance matrix could have an arbitrarily 
large eigenvalue ratio if the dimensions are chosen appropriately. Hence one could not 
examine the eigenvalue ratio in this case and conclude that the covariance matrix was 
diagonal. It is more appropriate to consider the eigenvalue ratio for the correlation 
matrix as an indicator of the condition of the matrix for inversion because in this case no 
dimensions are involved. This ratio is not readily available in the computer program 
used in this analysis, but the determinant of the correlation matrix is available and is 
used to indicate the condition of the matrix for inversion. Actually, the determinant con- 
tains some information about the eigenvalues since the product of the eigenvalues is equal 
to the determinant and the sum of the eigenvalues of the correlation matrix is equal to the 
order of the matrix. If q is the order and DET is the determinant of the correlation 
matrix, then it can be shown that 


DET 

4 



where and’ ^j^iin maximum and minimum values of the eigenvalues. This 

inequality gives a set of bounds on the eigenvalue ratio, however very little information is 
contained in the upper bound when the order of the determinant is large. 

After 1 orbit of tracking, the determinants of the correlation matrices in figure 3 
are extremely small and indicate nearly singular matrices. The smallness of these 
determinants is largely due to the short-data-arc problem mentioned previously. At the 
end of 10 orbits of tracking, the determinants have increased significantly but they are 
still so small that the associated normal matrices cannot be inverted with single precision 
or 8-decimal arithmetic. However, it appears that by using double -precision 16-decimal 
arithmetic, a first estimate of the 11 gravity coefficients considered herein can be obtained 
with just 1 orbit of tracking data. 


Comparison of Data Types 

As stated previously, two basic tracking data t 3 q>es are considered herein — range 
and range rate. It is of interest to compare the relative advantages of these two data 
types. One such comparison was given in the discussion of figure 2 where it was pointed 
out that the differences between the standard deviations in the coefficients obtained by 
using either data type alone or simultaneously are not significant. An equally important 
comparison is that between the corresponding correlation matrices. This comparison 
was made after each orbit up to 10 orbits and it was found that, after the third orbit, the 
correlation matrices for the two data t 3 q>es and their simultaneous use were similar. 
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That is, gravity coefficients which were highly correlated in one matrix were usually 
highly correlated in the other matrix, and, as expected, combining the two data types did 
not lead to any significant reduction in the correlations. A comparison of the determi- 
nants of the three correlation matrices was made after each orbit and it was found that, 
after 3 orbits of tracking, they were essentially the same. The matrices given in figure 3 
after 10 orbits of tracking are typical of these results. It was concluded for the values of 
Op and Op and the set of 11 gravity coefficients considered herein, that range and range- 
rate measurements are basically the same data types in the sense that they produce esti- 
mates of the gravity coefficients which are not significantly different and that they result in 
similar correlation matrices having essentially the same determinental values and that 
no significant reduction in correlations due to the combination of data types was noted. 

A similar conclusion was drawn in reference 6 in the case of estimating the lunar satel- 
lite state. 

Since it has been shown that range and range rate are similar data types, the results 
presented in the remainder of this paper are for the use of range rate only. The rela- 
tive advantage of range and range rate is affected by the values of the orbital elements, 
primarily a. For larger values of the semimajor axis, range has a slight advantage 
over range rate as far as the accuracy of estimating the coefficients is concerned. 

Effect of Semimajor Axis 

The variation of the standard deviation in the estimates of the gravity coefficients 
after 10 orbits of tracking with the semimajor axis are presented in figure 4. In order 
to obtain the standard deviations of the estimates of these 11 coefficients, normal mat- 
rices were formed and inverted after each orbit, up to 10 orbits, for each value of the 
semimajor axis investigated with the elements i, w, e, and tQ held constant. The 
semimajor axis was varied from 2 686 kilometers to 5 000 kilometers. Even though the 
period of the orbits changed as the semimajor axis chained, the number of range-rate 
measurements per orbit was held to 26. These measurements were assumed to be equally 
spaced in time over the orbital period and to have a standard deviation of 0.002 meter 
per second. Many of the curves in figure 4 have several local maximums and minimums 
which are due to changes in correlations between certain of the coefficients. In general, 
the standard deviation for a particular gravity coefficient increased when there was a 
significant increase in the correlation coefficients for that particular gravity and 
decreased with significant decreases in the correlations. There are several pairs of 
curves in figure 4 which are very similar, for instance those representing C 2 2 
S 2 !• This similarity is due to the high correlations between the two gravity coefficients. 
In this particular case, the two coefficients remained highly correlated over the entire 
range of semimajor axis. 
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With few exceptions the standard deviations presented in figure 4 increased as the 
semimajor axis increased. This trend is as expected since the partial derivatives which 
comprise the A matrix are inversely proportional to the n + c power of the semimajor 
axis where n is the degree of the coefficient and c is a constant. As the semimajor 
axis increased, the numerical values of the derivatives decreased and resulted in less 
sensitive derivatives; hence, each observation contains less information about the gravity 
coefficients. Because the information about the coefficients decreased, it is natural to 
expect a less accurate knowledge of the coefficients. Exceptions to the trend of increasing 
standard deviations can be seen in figure 4 and these exceptions are due to the previously 
noted changes in correlations between certain of the gravity coefficients. A good example 
of these exceptions is the curves representing C 2^1 and S 2 2 * Note that even these 
curves increase up to about 4 000 kilometers, drop off sharply between 4 000 kilometers 
and 4 750 kilometers and then begin to increase again. The sharp dropoff or steep neg- 
ative slope is due to a large reduction in the correlation between these two parameters. 

The slope of the curves in figure 4 would also be expected to increase as the degree 
of the coefficient increases because the partial derivatives in the A matrix are inversely 
proportional to the n + c power of the semimajor axis. This particular trend can be 
noted in figxire 4. 

Although a few correlations increased or decreased as the semimajor axis was 
varied over the given range, there were only a few very high correlations which changed 
appreciably. Also a comparison of the determinants of the correlation matrices after 
each orbit of tracking for each value of the semimajor axis investigated showed no signif- 
icant differences, and, hence, the condition of the normal matrix for inversion did not 
change appreciably over the entire range of semimajor axes. 

Effects of Inclination 

The inclination i of the satellite orbital plane with respect to the earth-moon 
plane has proven to be a significant parameter in the estimation of certain of the gravita- 
tional coefficients; namely, the even zonal coefficients. The reader is reminded that the 
assumption has been made that the earth-moon plane and the lunar equatorial plane coin- 
cide. The standard deviations in the estimates of the 11 coefficients after 10 orbits of 
tracking are plotted as functions of inclination in figure 5. Again the data presented in 
this figure were obtained by holding five of the orbital elements (a, co, e, and to) 
at the nominal values and var 3 dng the element of interest, in this case i. Note that the 
accuracy of estimating the two odd zonals C 3 q and Cg q is relatively independent of 
the inclination. The rapid decreases exhibited in some of the curves, especially the 
curves for the even zonal coefficients, are due to the reduction in correlations between 
certain of the coefficients with inclination. Likewise, increases in some of the curves 
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such as those for certain of the tesseral coefficients and for the even zonal Cq q 
due to increases in the correlations between certain coefficients with inclination. Although 
some of the tesserals may have correlations which are dependent on inclination, it is 
usually the correlation between pairs of even zonal coefficients which are most signifi- 
cantly influenced by the inclination. Examples of these correlations are shown in figure 6 
where the correlation matrices after 10 orbits of tracking are shown for inclination angles 
of 2°, 30°, and 60° and for a combination of 2° and 60°. Note that the correlations 
between certain pairs of even zonal coefficients such as C2^o ^4 o through 

zero at some inclination angle between 30° and 60° which is an indication that there exists 
some optimum angle between these limits for separation of these coefficients in the solu- 
tion vector. This result leads one to expect that the coefficients which have correlations 
that are highly dependent on inclination angle can be separated in the solution vector by 
combining tracking data from various satellites having different orbital inclinations. 

A comparison of the determinant of each of the correlation matrices presented in 
figure 6 shows that there also exists an optimum inclination angle, as far as the condition 
of the normal matrix is concerned. That is, there exists an inclination angle at which 
the determinant of the correlation matrix is a maximum. For the results presented, 
this angle is in the neighborhood of 40°. However, the important point is that the normal 
matrix might be better conditioned for inversion by combining data from various satellites 
having different orbital inclinations. 

A comparison of the determinants of the correlation matrices obtained by combining 
trackii^ data from satellites having different inclinations with those obtained by using 
tracking data from a single satellite is given in figure 7. Here it was assumed that 
26 range-rate measurements were made each orbit for 10 consecutive orbits of two dif- 
ferent satellites having orbital inclinations of 2° and 60°. The curve representing the 
results for the combined tracking data was obtained by adding the normal matrices for 
the two satellites after each orbit and inverting the resulting matrix to obtain the combined 
covariance matrix. For example, 10 orbits represented by the combined-data curve would 
include 5 orbits at i =* 2° and 5 orbits at i = 60°. It should be pointed out that the time 
interval of tracking in this case would be one-half that for 10 orbits of a single satellite. 
The results given in figure 7 are a clear indication that, in most cases, the normal matrix 
can be better conditioned for inversion by combining trackii^ data taken from satellites 
having different orbital inclinations. For instance, at the end of 10 orbits of tracking, 
the determinant of the correlation matrices for 2° or 60° is of the order 10“^® whereas 
the combined correlation matrix determinant is of the order 10“®. This represents a 
very significant increase in the determinant and hence a better conditioned normal matrix. 

A point of interest is that the determinant associated with the 60° inclination curve 
decreases with tracking time. This result was unexpected and is partly due to the 
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increasing correlations of the second-degree tesseral coefficients with certain of the 
other coefficients. The correlation matrix for the combined set of tracking data is given 
in figure 6 . By comparing this matrix with the others in figure 6 one can see that those 
correlations which have been pointed out as being dependent on inclination, the even zonals, 
are greatly reduced when the tracking data are combined. For example, the correlation 
between Cq q and C 2 q 2° and 60° inclination is 0.98 and -0.87, respectively, but 
with the combined tracking data, this correlation becomes 0.14. From these results it 
appears, as expected, that one way to reduce high correlations between certain of the 
gravitational coefficients and to obtain a better conditioned normal matrix would be to 
combine data from satellites having different orbital inclinations. 

Effect of Nodal Position 

One of the most important parameters in the determination of the gravitational 
coefficients is the node angle, which is the angle between the earth-moon line and the line 
of intersection of the satellite orbital plane with the earth-moon plane. The standard 
deviations in the estimates of the coefficients after 10 consecutive orbits of tracking at 
various nodal positions are presented in tabular form in figure 8 . These data were not 
plotted because the interval between the values of fi was not small enough to allow a 
smooth curve to be fitted to the points. It should be noted that the value of shown in 
figure 8 is the value at the beginning of the tracking interval and that at the end of 10 orbits 
Q has decreased approximately 19°. Of the 11 coefficients, the standard deviations for 
the two fourth-degree coefficients varied by approximately one order of magnitude where- 
as the remaining coefficients varied less than one order of magnitude. However, if the 
standard deviations are plotted as functions of the nodal positions, the resulting curves 
will contain a number of local maximums and minimums. This result is due to the strong 
dependence of the correlations between certain of the gravitational coefficients on the 
nodal positions. In regions of local maximums one would expect these gravitational coeffi- 
cients to be highly correlated, whereas in the regions of the local minimums the correla- 
tions have been reduced considerably. Examples of nodal-dependent correlations can be 
seen in figure 9 where the correlation matrices after 10 orbits of tracking are presented 
for two different nodal positions. For instance note that when equals 30° the corre- 
lation between C 2 j and S 2 2 ^•’ad the correlation between C 2 2 ^2 j are -0.97 

and 0.91, respectively, whereas for the case of equal 110 ° these two correlations 
have been reduced to 0.15 and 0.12, respectively. Although the correlations between these 
pairs of tesseral coefficients decreased considerably when was increased by 80°, 
there are correlations which increase considerably under the same circumstances, for 
example, the correlation between 82 ^ and 822 * Usually the nodal-dependent correla- 
tions are between pairs of tesseral coefficients. One implication of these nodal-dependent 
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correlations is that perhaps tracking should be spread throughout the lunar month to take 
advantage of viewing the lunar satellite orbit from all angles. 

In any real orbit determination process it is unlikely that a month of continuous 
tracking data from one particular satellite will be available. A more likely tracking 
schedule will be one which includes intermittent tracking throughout the month and, hence, 
tracking at various nodal positions. Since it has been shown in these analyses that the 
nodal position during the tracking interval is a very important parameter, it is of interest 
to compare the standard deviations of the gravity coefficients and the determinant of the 
correlation matrix for different tracking schedules. One such comparison is given in fig- 
ure 10 where the determinant of the correlation matrix is plotted as a function of the num- 
ber of orbits tracked for five different tracking schedules. In all five cases it was 
assumed that the tracking station would be available for tracking a lunar satellite for 20 
orbits dviring a period of 1 month. Curve (1) represents the results obtained by assuming 
the satellite had been tracked 1 orbit every 36 hours for a month while curve (2) is for a 
trackir^ schedule of 4 orbits every seventh day for a month, and curve (3) is for 2 orbits 
every third day for a month. For curves (4) and (5) it was assumed that the vehicle was 
tracked for 20 consecutive orbits with the nodal position at the beginnii^ of the tracking 
interval being 45° and 120°, respectively. The curves in figure 10 resulting from 20 con- 
secutive orbits of tracking at nodal positions of 45° and 120° represent approximately the 
upper and lower bounds for similar curves at all other nodal positions. The frequency of 
observations was 15 range-rate measurements per orbit. As far as the condition of the 
normal matrix is concerned, it is clear from figure 10 that tracking for 1 orbit every day 
and a half for a month is better than 4 orbits every 7 days or 2 orbits every 3 days for a 
month or for 20 consecutive orbits. However at the end of 20 orbits, tracking schedules 
(1), (2), and (3) result in nearly the same correlation matrix determinant which is approx- 
imately six orders of magnitude larger than the other two. This difference is significant 
and indicates that, by choosing the tracking schedule properly, the condition of the normal 
matrix for inversion can be improved. 

For tracking intervals less than 10 orbits, the standard deviations in the gravity 
coefficients associated with tracking schedule (1) are smaller than those associated with 
the other four tracking schedules. But again, as in the case of the determinant of the 
correlation matrix, the standard deviations after 20 orbits are not very different for 
tracking schedules (1), (2), and (3). However, the standard deviations associated with 
tracking schedules (4) and (5) are from one to three orders of magnitude larger than 
those of the other two schedules. Therefore, for a fixed allowable tracking time, the 
results indicate that the best tracking schedule would be one in which the observations are 
distributed over various nodal positions, as opposed to concentrating all the tracking at 
one nodal position. 
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Effect of Eccentricity 


The accuracy of estimating the gravitational coefficients was found to be a very 
strong function of the eccentricity. In figure 11 the standard deviations in the estimates 
of the 11 gravitational coefficients after 10 consecutive orbits of tracking have been 
plotted against the eccentricity. The same tracking schedule as previously noted was 
assumed, and again five of the orbital elements were fixed and the sixth, eccentricity, 
was varied over the given range. However, in this case, the nominal value of the semi- 
major axis was assumed to be 5 000 kilometers as opposed to 2 686 kilometers (which was 
used in all other cases). The large value of the semimajor axis was necessary to insure 
that, over the range of eccentricities investigated, the radius of periapse was never less 
than the radius of the moon. With the exception of the second-degree tesserals, the 
standard deviations in the estimates of the gravity coefficients varied by more than an 
order of magnitude and in the case of the fifth-degree zonal by three orders of magnitude 
over the given range of eccentricities. Some of the large variations in the standard devi- 
ations shown in figure 11 are due to reductions in the correlations between certain of the 
gravity coefficients. However most of them are due to the complicated way in which the 
eccentricity enters into the equations for the partial derivatives of the observations with 
respect to the gravity coefficients. For an eccentricity larger than approximately 0.2, 
all the gravity coefficients were better estimated more accurately as the eccentricity 
increased. A point of interest which was noted in the analysis of the eccentricity varia- 
tion was that, as far as the accuracy of estimating the gravity coefficients is concerned, 
range data had a slight advantage over range-rate data. Similar results were noted in 
reference 6 in the case of estimating the lunar satellite state, and they are due to the 
large value of the semimajor axis chosen for the nominal orbit. As was noted previously, 
at large semimajor axes range is a slightly better data type than range rate. 

Although it has been stated that certain correlations change with eccentricity, only 
a few high correlations have been noted to be dependent on eccentricity. In figure 12 the 
correlation matrices after 10 orbits of tracking are shown for an eccentricity of 0.01 
and 0.6. Several examples of eccentricity-dependent correlations are given in figure 12, 
for instance the correlation between C 2 q and q and the correlation between C 2 2 
and S 2 ^i- The correlations between these two pairs of gravity coefficients have also 
been pointed out as being dependent on inclination and nodal position, respectively. The 
correlations which were dependent on inclination were usually between pairs of zonal 
coefficients whereas the nodal-dependent correlations were between pairs of tesseral 
coefficients. No such consistency has been noted for the eccentricity-dependent correla- 
tions. Again, one would expect that certain of the gravity coefficients which are highly 
correlated could be separated in the solution vector by combining tracking data taken 
over satellite orbits havii^ different eccentricities. 
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Normal matrices were formed and inverted by using combined sets of assumed 
rai^e-rate measurements of lunar satellites having different orbital eccentricities. This 
was done to determine whether combining the range-rate measurements from two different 
eccentric orbits would reduce high correlations and improve the condition of the normal 
matrix for inversion. The normal matrices were formed and inverted after each orbit 
up to 10 orbits with 26 range-rate measurements per orbit assumed for each satellite. 
Hence, after 10 orbital periods for each satellite, the normal matrix would contain 20 
orbits of data. For the cases investigated, the condition of the normal matrix was not 
significantly improved for inversion, and no large reductions in correlations between 
highly correlated coefficients were noted. In fact some of the correlation coefficients 
increased when the tracking data were combined, and the determinant of the correlation 
matrix for the combined data sets was sometimes smaller than for the individual set. 

In several cases it was noticed that the standard deviations associated with the smallest 
determinant were smaller than those associated with a larger determinant. It appears 
that combining tracking data taken from orbits having different eccentricities is not a 
very effective method of reducing high correlations or of improving the normal matrix for 
inversion. 


CONCLUDING REMARKS 

In the real case of determining the lunar gravitational coefficients by using actual 
range and range-rate measurements of a lunar satellite, analyses of the types given in 
this report would be useful when choosing the set of gravity coefficients which can be best 
estimated with the available tracking data. The results presented herein should aid in 
pointing out some of the more important parameters for future research in the area of 
gravitational-parameter determination. 

The direct method of analysis for determining the gravitational coefficients of the 
spherical harmonics of the lunar potential function has been analyzed and seems to be a 
useful tool. A parametric study was made to determine the effects of the orbital elements 
on the accuracy of estimating the lunar gravitational coefficients and the condition of the 
normal matrix for inversion. The standard deviations in the estimates of certain of the 
gravity coefficients were very strong functions of the nodal position, eccentricity, inclina- 
tion, and the semimajor axis. Some of the high correlations associated with these esti- 
mates were found to be dependent either on the length of the tracking-data arc or on certain 
of the orbital elements. Usually the high correlations which were dependent on inclination 
were between pairs of even zonal coefficients whereas the nodal-dependent correlations 
were between pairs of tesseral coefficients. Several eccentricity-dependent correlations 
were noted but they were not associated with any particular set of coefficients. A com- 
parison of the results obtained by using range and range-rate measurements alone and by 
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using these measurements simultaneously was made to determine the relative advantages 
of the two data types. The correlation matrices resulting from either data type alone or 
from their simultaneous use were similar, and the determinant of the combined normal 
matrix was not significantly different from the individual normal matrices. Hence, it was 
concluded that range and range rate are basically the same data types. It appears that 
the best tracking schedule is one in which the allowable tracking time is distributed 
throughout the month. Finally, it was found that the condition of the normal matrix for 
inversion could be improved, and high correlations between the even zonal coefficients 
could be reduced when tracking data from lunar satellites having different orbital incli- 
nations are combined to estimate the gravity coefficients. 

It should be noted that the results presented herein are restricted to the particular 
set of 11 gravity coefficients investigated and if this set were altered, the results may be 
different. 

Langley Research Center, 

National Aeronautics and Space Administration, 

Langley Station, Hampton, Va., June 24, 1966. 
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Figure 3.- Correlation matrices after first and tenth orbit. 
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Figure 6.* Correlation matrices after 10 orbits for different inclinations. 
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Figure 8.- Standard deviations of estimates of lunar gravitational coefficients after 10 orbits of tracking at various nodal positions 

throughout month. 
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Figure 9.- Correlation matrix after 10 consecutive orbits for two different nodal positions. 
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Figure 12.- Correlation matrices for two different eccentricities. 
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